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The low-temperature states of bosonic fluids exhibit fundamental quantum effects 
at the macroscopic scale: the best-known examples are Bose-Einstein condensation 
(BEC) and superfluidity, which have been tested experimentally in a variety of differ- 
ent systems. When bosons are interacting, disorder can destroy condensation lead- 
ing to a so-called Bose glass. This phase has been very elusive to experiments due 
to the absence of any broken symmetry and of a flnite energy gap in the spectrum. 
Here we report the observation of a Bose glass of fleld-induced magnetic quasiparti- 
cles in a doped quantum magnet (Br-doped dichloro-tetrakis-thiourea-Nickel, DTN). 
The physics of DTN in a magnetic fleld is equivalent to that of a lattice gas of bosons 
in the grand-canonical ensemble; Br-doping introduces disorder in the hoppings and 
interaction strengths, leading to localization of the bosons into a Bose glass down to 
zero fleld, where it acquires the nature of an incompressible Mott glass. The transition 
from the Bose glass (corresponding to a gapless spin liquid) to the BEC (corresponding 
to a magnetically ordered phase) is marked by a novel, universal exponent governing 
the scaling on the critical temperature with the applied fleld, in excellent agreement 
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with theoretical predictions. Our study represents the first, quantitative account of the 
universal features of disordered bosons in the grand-canonical ensemble. 

PACS numbers: 03.75.Lm, 71.23.Ft, 68.65.Cd, 72.15.Rn 

Introduction, Disorder can have a very strong impact on quantum fluids. Due to their v^ave-like 
nature, quantum particles are subject to destructive interference w^hen scattering against disordered 
potentials. This leads to their quantum localization (or Anderson localization), vv^hich prevents e,g. 
electrons from conducting electrical currents in strongly disordered metals [IJ, and non-interacting 
bosons from condensing into a zero-momentum state ]Q^. Yet interacting bosons represent a matter 
w^ave v^ith arbitrarily strong non-linearity, w^hose localization properties in a random environment 
cannot be deduced from the standard theory of Anderson localization. For strongly interacting 
bosons it is knovv^n that Anderson localization manifests itself in the Bose glass: in this phase the 
collective modes of the system - and not the individual particles - are Anderson-localized over 
arbitrarily large regions, leading to a gapless energy spectrum, and a finite compressibility of the 
fluid Ell]]. Moreover nonlinear bosonic matter weaves can undergo a localization-delocalization 
quantum phase transition in any spatial dimension w^hen the interaction strength is varied [[3ll3|; the 
transition brings the system from a non-interacting Anderson insulator to an interacting superfluid 
condensate, or from a superfluid to a Bose glass. Such a transition is relevant for a large variety 
of physical systems, including superfluid helium in porous media [6J, Cooper pairs in disordered 
superconductors [7J, and cold atoms in random optical potentials ElU. Despite the long activity 
on the subject, a quantitative understanding of the phase diagram of disordered and interacting 
bosons based on experiments is still lacking. 

Recent experiments have demonstrated the capability of realizing and controlling novel Bose 
fluids made of quasiparticles in condensed matter systems [[9l[l0l. In this context, a prominent 
place is occupied by the equilibrium Bose fluid realized in quantum magnets subject to a mag- 
netic field [lOJ, in v^hich disorder can be introduced in a controlled w^ay by chemical doping, 
leading to novel bosonic phases [|TT1 [12]|. The ground state of such systems v^ithout disorder 
and in zero field corresponds to a gapped bosonic Mott insulator. Extra bosons can be injected 
into the system by applying a critical magnetic field that overcomes the gap, and that drives a 
transition to a superfluid state (magnetic Bose-Einstein condensate - BEC). Such a state corre- 
sponds to an XY antiferromagnetic state of the spin components transverse to the field. Here 
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we investigate the Bose fluid of magnetic quasiparticles realized in the model S = 1 compound 
NiCl2-4SC(NH2)2 (dichloro-tetrakis-thiourea-Nickel, DTN) [13J via experiments (AC magnetic 
susceptibility, DC magnetization and specific heat), and large-scale quantum Monte Carlo (QMC) 
simulations. Disorder is introduced by Cl^Br substitution, which, as we will see, leads to ran- 
domness in the bosonic hoppings and interactions. We observe a Bose glass in two extended 
regions of the temperature-magnetic field phase diagram of Br-doped DTN. The gapless nature 
of the Bose glass manifests itself in a finite uniform magnetic susceptibility (corresponding to the 
compressibility of the quasiparticles), and in a non-exponential decay of the specific heat at low 
temperature, probing the low-energy density of states. This state extends down to zero field: in this 
limit the compressibility/susceptibility vanishes, while the spectrum remains gapless, giving rise 
to a Mott glass. We investigate the thermodynamic signatures of the Mott and Bose glasses, and 
the Bose-glass-to-superfluid transition, characterized by a novel universal exponent for the scaling 
of the condensation temperature with applied field. 

Magnetic properties of pure DTN. The magnetic properties of pure DTN are those of antiferro- 
magnetic S = 1 chains of Ni^+ ions, oriented along the crystallographic c-axis, and coupled trans- 
versely in the a6-plane lfT3] - [T5ll . A strong single-ion anisotropy D is present, with an anisotropy 
axis z corresponding to the c-axis. The magnetic Hamiltonian reads 



where Jc = 2.2 K is the antiferromagnetic coupling for bonds {ij)c along the c-axis. Jab = 0.18 
K is the coupling for bonds {lm)ah in the ah plane, and = 8.9 K is the single-ion anisotropy. 

= 2.26 is the gyromagnetic factor g along the c-axis. In zero field, the large single-ion anisotropy 
D forces the system into a quantum paramagnetic state with each spin close to its \ms = 0) state. 
Mapping the aS = 1 spin states onto bosonic states with occupation n = ms + 1, the quantum 
paramagnet corresponds to a Mott insulator of bosons with n = 1 particles per Ni site, and with 
a gap /\ D — 2Jc — 4: Jab + 0{Jl/D) for the addition of an extra boson. A magnetic field 
exceeding the value = /S./ gjiB ^ 2. IT is able to close the spin gap and to create a finite 
density of excess bosons that condense into a magnetic BEC (see Fig. [if a)). The appearance of 
excess bosons translates into a finite magnetization along the field axis; their long-range phase 
coherence translates into long-range XY antiferromagnetic order transverse to the field. Long- 
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range order persists up to a critical condensation temperature which, for H > H^i\ scales 
with the applied field as ^ \H — Here = 2/3, as predicted by mean-field theory 

for a diluted gas of excess bosons, and as measured with very high accuracy down to ImK lfT6ll . 
When the magnetic field is increased further, the spins are brought to saturation for H = H^^^ = 
{D + 4 Jc + 8Jab)/gf^B = 12.6 T, and the system transitions from a BEC to a perfect Mott insulator 
with n = 2 particles per site. Correspondingly, the BEC critical temperature vanishes as Tc ^ 

Phase diagram of Br-doped DTN. We have measured the critical temperatures and fields for 
magnetic BEC in Ni(Cli_^Br^)2-4SC(NH2)2 (Br-DTN) with x = 0.08 ± 0.005 doping by measur- 
ing the AC susceptibility at low frequencies and the specific heat ifTTll . AC and DC susceptibility 
measurements are performed at fixed temperature and varying fields, and they show a step-like 
increase/decrease corresponding to the critical field for BEC, similar to the pure sample lfT6l [see 
Fig. [3]^a)-(b)]. The main difference compared to pure DTN is that - at low temperatures - the 
upper and lower edge of the steps are rounded by disorder; as we will see below, this rounding is 
a fundamental indication of the nature of the phases connected by the transition. An independent 
estimate of the critical BEC temperature as a function of the field is obtained by the location of a 
sharp A-peak in the specific heat [Fig. [3]^c)]. The sharpness of the features corresponding to the 
BEC transition are quite remarkable, given the strong doping introduced in the system. Moreover, 
for temperatures below the A-peak the specific heat clearly follows a behavior, consistent with 
long-range XY antiferromagnetic order in 3D. Fig. [2] summarizes the experimental phase diagram 
of Br-DTN. Br-doping has a profound impact on the phase diagram of DTN: in particular both 
the lower and upper critical fields for the onset of magnetic BEC at T ^ are found to shift to 
lower values, Hd = 1.07(1) T and Hc2 = 12.16(1) T, as shown in Fig.[2| But most importantly 
the magnetic behavior of Br-DTN outside the BEC region is completely different compared to the 
pure system. In the pure system, the ground state outside the magnetic BEC is a Mott insulator 
with a large spin gap A away from the critical fields. This leads to an exponential suppression of 
the specific heat at low temperatures fc^T < A as Cy ^ exp [-A/{kBT)], as shown in Fig.jsfd), 
and to a similarly vanishing susceptibility for T ^ 0. On the contrary, for x = 0.08, we ob- 
serve that the susceptibility is finite for H > Hc2, and it even exhibits a strong satellite peak for 
H ^ 13.5 T. The susceptibility vanishes only for H = Hg ^ 17 T, corresponding to the satura- 
tion field of the entire sample, which is pushed to a much higher value than in the pure sample 
(where Hg = H^2^ = 12.6 T). In the region H < Hd we observe that the specific heat exhibits 
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FIG. 1: Sketch of the bosonic phases of DTN and Br-doped DTN: In the undoped case, an increasing 
magnetic field along the c-axis drives the system from a Mott insulating (MI) phase (a) to a BEC phase (b) 
by injecting delocalized excess bosons (indicated in cyan) on top of the MI background at density n = 1; in 
the doped case, an arbitrarily weak magnetic field can inject extra bosons in the Br-rich regions (indicated 
by the orange bonds) which are localized and incoherent in the (low-field) Bose glass phase (c) - their 
localized wavefunction is sketched by the light-blue lines; further increasing the magnetic field leads to the 
percolation of phase coherence via coherent tunneling of the excess bosons between the localized regions, 
giving rise to an inhomogeneous BEC (d); for strong magnetic fields H < Hc2 the spins away from the 
Br-bonds are close to saturation/double occupancy (represented in dark blue), and unpolarized spins / singly 
occupied sites, corresponding to bosonic holes, only survive in the Br-rich regions (e). These holes localized 
into disconnected, mutually incoherent states when entering the high-field Bose glass (f). 



a non-exponential decay, down to zero field [Fig.[3]^d)]. Therefore we can conclude that the non- 
magnetic phases for < H < Hg correspond to gapless bosonic insulators, which, as we will see, 
can be identified with a compressible Bose glass (for H > 0) and an incompressible Mott glass 
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FIG. 2: Experimental phase diagram of Br-doped DTN from specific heat and susceptometry, compared 
to QMC data; the following phases are represented: Bose-Einstein condensation (BEC), Bose glass (BG), 
and Mott glass (MG). The lilac regions represent the magnitude of the spin gap in the Mott insulating (MI) 
phase; (b) Experimental phase diagram of pure DTN (based on specific heat and the magnetocaloric effect, 
see Ref. [14J). 



(for H = 0). 

Modeling Br doping. Br-DTN can be successfully modeled theoretically by considering that 
Br substitution for CI affects the super-exchange paths associated with the Jc couplings, and it 
also distorts the lattice locally due to the larger atomic radius of Br with respect to CI. The 
disappearance of the spin gap down to H = and the upward shift of the saturation field suggests 
that Br doping locally strengthens the magnetic coupling Jc and lowers the anisotropy D. For 
simplicity we only consider that Ni-Cl-Cl-Ni bonds in DTN can be turned into Ni-Cl-Br-Ni or 
Ni-Br-Cl-Ni, and we neglect Ni-Br-Br-Ni bonds that represent only 0.6% of the total bonds for 
X = 0.08. We assign a value to the magnetic exchange coupling of the Br-doped bonds, and a 
value to the single-ion anisotropics of the Ni ion adjacent to the Br dopant. Note that for a doping 
concentration x, we have a fraction of 2x doped bonds, given that each bond can accommodate a 
Br dopant on two different CI sites. We then use and as fitting parameters of the full low- 
temperature magnetization curve in Fig. [3]^a), which is calculated using QMC simulations ifTTll . 
We find an extremely good agreement between experimental data and simulation for ^ 2.35 Jc 
and ^ D/2, giving us confidence that we are able to quantitatively model the fundamental 
microscopic effects of doping in Br-DTN. Indeed the critical temperature for BEC, extracted from 
a finite-size scaling analysis of the simulation data with doping x = 0.075 II17I1 . is in remarkable 



7 



1 

0.8 



pure DTN 
— Br-DTN / 
QMC (doped case) /; 

// 

^''^ 0.1 




5 10 \5 
15 20 




- T=l mK 

- 70 mK 

- 150 mK 
200 mK 
350 mK 

- 500 mK 

- 787 mK 



- T=l mK 

- 100 mK 

- 150 mK 
250 mK 
400 mK 

- 600 mK 



0.4 0.6 
T (K) 




0.5 1 1.5 2 
H(T) 



12 13 14 

H(T) 



10 

r^lO 
^10 

|io- 



^^10' 

uio' 

10"^ 
10'" 
10'^ 



' ' ' ' 'QXp.(x=0) ' ' : 

- o exp. (x=0.08j^ 

- ■ theory (LGM)i 

^ ^ 

r H=OT : 

- 1 1 1 1 1 1^ 1 1 1 1 > 1 1 : 


r H=0.5T ^^p 

- 1 1 1 1 1 1 1 1 1 1 1 1 1 1 - 


r//=0.75T 1 

~ 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 ' 


fV-J 

! \ ! 
: \ : 

r H=IT ^ 1 



r ' (k') 



FIG. 3: (a) Magnetization curve of Br-DTN at T = 19 mK, compared to QMC results, and to pure DTN 
magnetization (measured at T = 16 mK). In the inset we show the DC susceptibility curve, obtained by 
differentiating the magnetization; (b) AC susceptibility of Br-DTN at frequency / = 88.7 Hz close to the 
lower and upper critical fields. The curves have been shifted with respect to one another for readability 
purposes. The arrows indicate the appearance of sharp kinks at higher temperatures; (c) Specific heat of Br- 
DTN from i7 = T to i7 = 2 T; (d) Specific heat of Br-DTN in the MG and BG phase for H < Hd ^ I 
T, showing a non-exponential decay as T ^ 0; a comparison is made to the theory predictions based 
on the local-gap model (LGM), and to the data for pure DTN; in the upper-left and lower-right panels, 
the blue dashed line is a fit of the pure-DTN data to Aexp{—A{H)/kBT) where A is a constant and 
A{H) /kB = 9f^B{H^^ - H) /ks = 3.16 K for = and 1.64 K for = 1 T 



quantitative agreement with the experiment, as shown in Fig. |2]^a). The critical fields estimated 
from simulations are = 1.199(5) T and = 12.302(5) T, slightly larger (by - 0.13 - 0.14 
T) than the experimental values. However the downward shifts of Hd and Hc2 with respect to the 
pure system are correctly captured. 

Bose and Mott glass. The anisotropy D acts as a repulsion term inhibiting two bosons from 
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sitting on the same site, while the couphng Jc controls the kinetic energy of the bosons. In partic- 
ular we find numerically ifTTl that a model in which all c-axis bonds contain a Br dopant (leading 
to a couplings and to an anisotropy D' on one of the two sites connected by the bond), is in a 
BEC phase (XY ordered phase) even in zero field. This means that the Br-rich regions in DTN, 
characterized by the Hamiltonian parameters and D' , behave locally as mini-BECs, and they 
are locally gapless. Strictly speaking, Br-rich regions will have a residual gap due to their finite 
size. However, the statistical distribution of sizes has no upper bound, so that the distribution of 
local gaps has no lower bound, and consequently Br-DTN is globally gapless even in zero field. 
The corresponding bosonic phase is therefore a gapless insulator with spin inversion symmetry 
along field axis, and resulting in a commensurate boson density n = 1. This represents to our 
knowledge the first experimental realization of a Mott glass [[HI [HI, which has the peculiar as- 
pect of being incompressible (namely of having a vanishing susceptibility at T = 0) despite being 
gapless [El [m. As soon as a field is applied to the system, it can immediately inject excess 
bosons, which localize a la Anderson in the Br-rich regions, resulting in a paradigmatic example 
of a Bose glass (Fig.[T]^c)). In the spin language, spins in the Br-rich regions acquire a finite mag- 
netization along the field and their transverse components correlate antiferromagnetically over a 
finite range, but the local phase of the antiferromagnetic order is different from region to region 
so that the system remains globally paramagnetic. Long-range phase coherence of the local order 
parameters - corresponding to the local phases of the bosonic wavefunction - is established only 
when the localized states of the bosons grow enough under the action of the applied field as to 
overlap, leading to coherent tunneling of bosons between neighboring localized states (Fig.[TJd)). 
The resulting phase is a highly inhomogeneous BEC [[20l. 

We can quantitatively test the picture of bosons localized in Br-rich regions against the thermo- 
dynamic behavior of Br-DTN by using a simplified local- gap model (LGM). Within this model 
|[T7]| . the low-temperature and low-field behavior of the system is reduced to that of a collection of 
three-level systems, corresponding to a local longitudinal magnetization m^-^tot = 0, ±1 for each 
localized state. There is a finite-size gap Aat ^ c/N (for zero field) between the ms-^tot = ground 
state and the m^' tot = ±1 excited states, where is the number of sites in the Br-rich cluster. The 
low-temperature specific heat in zero field can then be predicted analytically to be 



where t = ksT / Jc and = log(2x); this expression displays the stretched exponential behavior. 




(2) 
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that uniquely characterizes the Mott glass [fT9l . The c parameter, and an overall prefactor, are 
used as fitting parameters of the experimental data in zero field, leading to an extremely good fit, 
as shown in Figs. |3]^d) and|4} Notably, no further adjustable parameters are necessary to fit the 
finite-field data, displayed in Fig. |3]^d), which also show a remarkable agreement with the theory 
prediction upio H ^ Hd. 




FIG. 4: Specific heat in zero field, displaying the characteristic Mott glass scaling. The solid lines are 
theoretical predictions based on the LGM, with parameter c = 3.02. 

For H Hc2 the magnetization approaches the value m^, = 1 — 2x 0.84, where all spins not 
connected to a Br-doped bond are polarized - and indeed Hc2 lies very close to the polarization 
field H^2^ of pure DTN. The full polarization of the Br-poor regions leads to a pseudo-plaitm in 
the magnetization aim^rrix (pseudo because it still exhibits a small finite slope). This feature 
corresponds to the high-field BG phase, which is characterized by the localization of bosonic holes, 
or singly occupied sites with ms = 0, in a background of doubly occupied sites with ms = 1 
(Fig. [T]^f)). Such holes persist up to the saturation field Hg, which is the field necessary to fully 
polarize a homogeneous system with couplings and anisotropics everywhere. The step-like 
feature in the magnetization at the upper bound of the pseudo-plateau is therefore induced by the 
saturation of the Br-rich clusters, and it is smeared due to the fact that such clusters have random 
geometries and therefore a distribution of local saturation fields, upper bounded by Hg. 

Thermal percolation crossover The physics described so far is valid only for low temperatures. 
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FIG. 5: Scaling of the critical temperature with the distance from the T = critical fields, exhibiting a 
crossover between various exponents. The dashed and dotted lines indicate a fit to the form a|i7 — 77^1(2) 
and a\H — i?ci(2)l^^^ respectively, while the solid line is a fit to a'\H — i?ci(2)l'^' with the resulting 
exponent indicated in the figure (a and a' are fitting parameters). The leftmost panels show the critical 
line extracted from the AC susceptibility, and the rightmost ones the critical line obtained from the QMC 
simulations. 

As the temperature is increased above ^ 200 mK, the bosons that were localized in the Bose Glass 
state thermally delocalize and proliferate; this leads to a thermal percolation of their density pro- 
file (corresponding to the longitudinal magnetization profile) throughout the sample [l20]| . Thus a 
more ordinary paramagnet forms at higher temperatures and the nature of field-driven transition 
into the BEG phase also changes fundamentally. Indeed at temperatures below the thermal perco- 
lation crossover, the BEG transition occurs as sketched in Figs.[T]^c)-(d) and Figs.[T]^e)-(f), namely 
by coherent tunneling of bosons between localized states, resulting in a highly inhomogeneous 
BEG phase. This picture changes above the thermal percolation crossover. Now in the normal 
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phase, the bosons move incoherently on a pre-percolated network of magnetized sites, and their 
BEC transition upon increasing the field corresponds therefore to condensation on a random 3D 
percolated lattice, which is fully analogous to condensation on a regular 3D lattice. Signatures of 
the thermal percolation crossover can be found in the critical behavior of the AC susceptibility: At 
low temperatures (T < 200 mK) it exhibits a rounded shoulder for H > Hd and H < Hc2, and at 
higher temperatures it shows a sharp kink - analogous to what is observed in the pure system lfT6ll 
(seeFig.[3]:b)). 

But the most dramatic signature of the thermal percolation crossover is observed in the scaling of 
the critical temperature with the applied field, shown in Fig. [5} Plotting vs. \H — i/ci(2)| on a 
log-log scale, we clearly observe a kink separating two different scaling regimes. At high temper- 
atures (T > 200 — 300 mK) the field-dependence of Tc is essentially consistent with a pure-system 
scaling for low temperatures, Tc ^ \H — i/ci(2)|'^ with (p = 2/3, or with a pure-system scaling 
for intermediate temperatures with = 1/2, as observed in other magnetic BEC systems [21 J . 
At low temperatures, the scaling exponent crosses over to novel values, = 1.1(2) (close to Hd) 
and (/) = 1.1(1) (close to i/c2), which are consistent within the error (see [fTTl for a discussion of 
the estimate of cj)). Moreover, these novel scaling exponents are consistent as well with the values 
extracted from our QMC simulations (cp = 1.06(9) and 1.2(1) close to Hd and Hc2 respectively). 
Simulations also show a rough quantitative agreement for the crossover temperature range. Most 
remarkably, a consistent value of the exponent (/) at low temperature is also observed theoretically 
for the magnetic Hamiltonian of DTN subject to a different type of disorder, namely site dilution 
II20II . We can therefore conclude that the low-temperature scaling of Tc exhibits a novel exponent 
0^1 — 1.1 which is a universal feature of the BG-BEC transition. 

Conclusions. We have performed a comprehensive experimental and theoretical study of the 
disordered and strongly interacting Bose fluid realized in a doped quantum magnet (Br-DTN) 
under application of a magnetic field. We provide substantial evidence of the existence of gapless 
insulating phases of the bosons - the Mott glass and the Bose glass - and we investigate for the first 
time the quantitative features associated with their thermodynamic behavior. These phases can be 
quantitatively described as a Bose fluid fragmented over an extensive number of localized states 
with variable local gaps, dominating the response of the system. The presence of a Bose glass leads 
to a novel and seemingly universal exponent governing the scaling of the critical temperature for 
the transition from Bose glass to BEC. The remarkable agreement between theory and experiment 
shows that Br-DTN is an extremely well controlled realization of a disordered Bose fluid, which 
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allows a detailed experimental study of the thermal phase diagram of disordered bosons in the 
grand-canonical ensemble. 

Acknowledgements. Work at the High Magnetic Field Laboratory at the Physics Institute of 
the University of Sao Paulo were supported in part by the Brazilian agencies FAPESP and CNPq. 
Measurements at the NHMFL High B/T and pulsed field facilities were supported by NSF Grant 
DMR 06541 18, by the State of Florida, and the DOE. Work at LANL was supported by the NSF, 
and the DOE's Laboratory Directed Research and Development program under 20100043DR. 
The numerical simulations have been performed on the computer facilities of the NCCS at the 
Oak Ridge National Laboratories, and supported by the INCITE Award MAT013 of the Office of 
Science - DOE. 



[1] B. Kramer and A. MacKinnon, Rep. Prog. Phys. 56, 1469 (1993). 

[2] L. Fallani, C. Fort, and M. Inguscio, Adv. At. Mol. Opt. Phys. 56, 1 19 (2008). 

[3] T. Giamarchi and H. J. Schulz, Phys. Rev B 37, 325 (1988). 

[4] M. P A. Fisher et ai, Phys. Rev. B 40, 546 (1989). 

[5] P B. Weichman, Mod. Phys. Lett. B 22, 2623 (2008). 

[6] P A. Crowell et al, Phys. Rev. B 55, 12620 (1997). 

[7] B. Sacepe et al, Nature Phys. 7, 239 (201 1). 

[8] L. Sanchez-Palencia and M. Lewenstein, Nature Phys. 6, 87 (2010). 

[9] H. Deng, H. Haug and Y. Yamamoto, Rev. Mod. Phys. 82, 1489 (2010), and references therein. 
[10] T. Giamarchi et al. Nature Phys. 4, 198 (2008), and references therein. 

[11] O. Nohadani et al, Phys. Rev. Lett. 95, 227201 (2005); T. Roscilde and S. Haas, Phys. Rev. Lett. 95, 

207206 (2005); T. Roscilde, Phys. Rev. B 74, 144418 (2006). 
[12] H. Manaka, A. V. Kolomiets, and T. Goto, Phys. Rev Lett. 101 077204 (2008); T. Hong, A. Zheludev, 

H. Manaka, and L.-P Regnault, Phys. Rev. B 81 060410 (2010). 
[13] V. S. Zapf et al, Phys. Rev. Lett. 96, 077204 (2006). 
[14] S. A. Zvyagin et al, Phys. Rev. Lett. 98, 047205 (2007). 

[15] The structure of DTN is actually that of two interpenetrating tetragonal lattices, which can be consid- 
ered effectively as decoupled, see S. A. Zvyagin et al, Phys. Rev. B 77, 092413 (2008). 
[16] L. Yin et al, Phys. Rev. Lett. 101, 187205 (2008); J. Low. Temp. Phys. 158, 710 (2010). 



13 

[17] Supplementary Material. 

[18] N. Prokof'ev and B. Svistunov, Phys. Rev. Lett. 92, 015703 (2004). 
[19] T. Roscilde and S. Haas, Phys. Rev. Lett. 99, 047205 (2007). 
[20] R. Yu et al, Europhys. Lett. 89, 10009 (2010). 
[21] N. Kawashima, J. Phys. Soc. Jpn. 73, 3219 (2004). 



SUPPORTING ONLINE MATERIAL 
Bose glass and Mott glass of quasiparticles in a doped quantum magnet 



Rong Yu/ Liang Yin,^ Neil S. Sullivan,^ J. S. Xia,^ Chao Huan,^ Armando Paduan-Filho,^ 

Nei F. Oliveira Jr.,^ Stephan Haas,"^ Alexander Steppke,^ Corneliu F. Miclea,^ Franziska 
Weickert,^ Roman Movshovich,^ Eun-Deok Mun,^ Vivien S. Zapf,^ and Tommaso Roscilde^ 

^Department of Physics & Astronomy, Rice University, Houston, TX 77005, USA 
^ Department of Physics and National High Magnetic Field Laboratory, University of Florida, Gainesville, FL 32611, USA 
^Instituto de Fisica, Universidade de Sao Paulo, 05315-970 Sao Paulo, Brasil 
Department of Physics and Astronomy, University of Southern California, Los Angeles, CA 90089-0484, USA 
^ Max-Planck Institute for Chemical Physics of Solids, Nothnitzer Str 40, 01187 Dresden, Germany 

^Condensed Matter and Magnet Science, Los Alamos National Lab, Los Alamos, NM 87545 
^ Laboratoire de Physique, Ecole Normale Superieure de Lyon, 46 Allee d' Italic, 69007 Lyon, France 



PACS numbers: 03.75.Lm, 71.23.Ft, 68.65.Cd, 72.15.Rn 




FIG. 1: (a) NiCb lattice in DTN; (b) effect of Br doping on the 
super-exchange coupling between Ni^"^ ions. 



SAMPLE PREPARATION AND EXPERIMENTAL 
METHODS 

DTN crystals were prepared by evaporation of an aqueous 
solutions of thiourea and nickel chloride. A slight excess of 
hydrochloric acid facilitates the growth of the samples. The 
solution for the doped crystals were prepared by thiourea, 
nickel chloride and hydrobromic acid. The solution temper- 
ature was kept around 35 degrees centigrade. The amount of 
Br/Cl atoms in the sample were determined by potentiometric 
analysis. The structure of the Ni-Cl lattice in DTN is shown 
in Fig. 1(a). The Ni atoms form two interpenetrating tetrago- 
nal systems which are effectively decoupled from a magnetic 
point of view [1]. Br substitution alters the bonds connecting 
two adjacent Ni atoms, and it is found to alter the magnetic 
coupling between the Ni ions (from Jc to J'^) as well as the 
on-site single-ion anisotropics (from D to D'). 

All magnetic measurements were made with the magnetic 
field applied along the tetragonal axis (c-axis) of the sam- 
ple. The AC susceptibility measurements were carried out 
using a PrNis nuclear refrigerator (down to 1 mK) and a 
15 T magnet at the High B/T facility of National High Mag- 
netic Field Laboratory (NHMFL). The experimental methods 



have been described elsewhere [2]. The AC susceptibility 
measurements allow to extract the critical field for magnetic 
BEC at fixed temperature, H = Hc{T), which provides the 
magnetic-field dependence of the critical temperature by in- 
version, T — Tc{H). The lower/upper critical field He is 
determined by a peak/dip in the first derivative of the AC sus- 
ceptibility with respect to the applied field. 

The DC magnetization was measured at the Physics Insti- 
tute of the University of Sao Paulo by using a plastic dilution 
refrigerator at 19mK. The main magnetic field was generated 
by a NbaSn superconducting magnet, with a maximum field 
of 17 T. The field sweep rate, [dB/dt], was 0.43 T/min. The 
magnetization was measured by the force method: the mag- 
netic force acting on the sample was produced by a magnetic- 
field gradient, generated by a superconducting gradient coil 
superimposed on the main magnetic field, and it was mea- 
sured by a capacitance technique. Specific heat was measured 
in a Quantum Design ^^He/^He dilution refrigerator down to 
50 mK using the thermal relaxation method. 



THEORETICAL METHODS 



Spin-to-boson mapping 



The model for the magnetic behavior of DTN and Br-doped 
DTN, Eq. (1) of the main text, can be exactly cast in the 
form of a Hamiltonian for interacting bosons via the Holstein- 
Primakoff spin-boson transformation for 5 = 1 spins: S~ = 
V2 y/l-n/2 b, = V2 6t^l -ri/2, and = n - 1. 
Here n = b'^b, and b, 6^ are triplet-boson operators satisfy- 
ing the commutation relation [6, 6^] = 1 and the constraint 
n = {0, 1, 2}. The magnetic Hamiltonian modeling Br-doped 
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DTN takes then the bosonic form 



n 



+ ^ J,,(n, - l){nj - 1) + ^ A(riz - 1)^ 
—qIIbH 77,^ + const. (1) 



Here J^^ = Jab for nearest-neighbor bonds {ij) lying in the 
a6 plane; for bonds along the c axis, Jij can take values Jc 
or = 2.3 Jc with probabilities 1 — 2x and 2x respectively, 
where x is the Br concentration. The single-ion anisotropy Di 
takes the value D' when the site is adjacent to a Br dopant, 
and value D otherwise (with D' = D/2). 

It is evident that Eq. (1) has the form of an extended Bose- 
Hubbard Hamiltonian describing nearest neighbor hopping, 
on-site and nearest neighbor interactions of bosons on a lat- 
tice. In particular, the magnetic couplings play the role of an 
hopping term for the bosons as well as of nearest neighbor re- 
pulsion. In Eq. (1) we have exploited the bipartite nature of 
the magnetic lattice to change the sign of the hopping term via 
a gauge transformation hi — 6^ on a sublattice. Compared 
to a standard Bose-Hubbard Hamiltonian, extra ^1 — ^ fac- 
tors appear in the hopping term, which provide the correct val- 
ues of the matrix elements of the spin operators and S~, 
and in particular they constrain the bosonic occupation num- 
ber. They do not alter the U(l) symmetry of the Hamiltonian, 
which is spontaneously broken by the BEC state (XY ordered 
state). Therefore Eq. (1) shares the same universal critical fea- 
tures of a standard Bose-Hubbard Hamiltonian. The single- 
ion anisotropy plays the role of an on-site repulsion which 
is particle-hole symmetric around the ri = 1 state, namely 
it penalizes equally the double occupancy n = 2 and the 
vacuum n = 0. The nearest-neighbor repulsion due to Jij 
is also particle-hole symmetric, and it favors a doubly occu- 
pied site neighboring an empty one. However the dominant 
D anisotropy makes this configuration very unlikely, so that 
the n.n. repulsion plays a marginal role. Br-DTN can thus be 
modeled using an extended Bose-Hubbard Hamiltonian with 
random hoppings, random on-site interactions, and random 
nearest-neighbor interactions, which are bimodally distributed 
in a spatially correlated way. 



been estimated via a careful finite- size scaling analysis of the 
spin correlation length, which is calculated via the disorder- 
averaged second-moment estimator [4]: 
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Here L is the linear dimension of the finite system, Xp refers 
to the unit vector along the p (p — a{b), c) crystal axis, and 
the ordering wavevector Q = (7r,7r,7r). [. . .]av denotes the 
disorder average. S{Q) is the static spin structure factor for 
the transverse spin components. 
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Near the critical point, the correlation length follows the scal- 
ing ansatz 



(4) 



which is verified by the QMC data shown in Fig. 2. At each 
temperature, we determine the critical field /ici(2) as the cross- 
ing point of the .^/L curves for different lattice sizes. 
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FIG. 2: Finite-size scaling of the correlation length at T ■ 
in the Br-doped DTN model with x = 0.075. 



172 mK 



Quantum Monte Carlo simulations 

To obtain the theoretical phase diagram of the model for Br 
doped DTN, we perform numerically exact quantum Monte 
Carlo (QMC) simulations based on the stochastic series ex- 
pansion algorithm [3]. The simulations are done on cubic 
lattices with up to IS"^ quantum spins and the simulation re- 
sults are averaged over 300 disorder realizations. For some 
selected temperatures and fields we have run simulations on 
lattice with 20"^ spins and averaging over 1200 disorder re- 
alizations. The critical temperatures and critical fields have 



We also applied the QMC simulations to study the phase 
diagram of the system with homogeneous couplings on all 
c-axis bonds and having anisotropics D' on one of the two 
sites connected by the bond. This models the idealized sit- 
uation (denoted as Br-only DTN) in which all all the bonds 
along the c-axis are of the Ni-Br-Cl-Ni type, and in which the 
position of the Br atom is randomly chosen between the two 
CI sites; therefore the anisotropics are randomized, taking on 
average the value D' on 3/4 of the sites and D on the reman- 
ing 1/4. Although not realizable in the experiments as such. 
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this model extends the behavior of the Br-rich regions in real 
Br-doped DTN to the thermodynamic limit, and therefore it 
unveils the local response of the Br-rich regions. Following 
a similar analysis of the correlation length as that described 
above, we find that Br-only DTN displays XY antiferromag- 
netic order up to a finite critical temperature even in zero mag- 
netic field, and it therefore exhibits a gapless spectrum. As 
shown in Fig. 3, the finite-size scaling of ^/L is consistent 
with a divergence in the thermodynamic limit at = and a 
finite temperature T = 68.8 mK, indicating that the system is 
long-range ordered at this temperature. We notice that this is a 
non-trivial result, given that we are dealing with a random- L) 
system. On the other hand we have studied separately the non- 
random situation in which all anisotropics take the value D, 
finding that even in this case the system displays long-range 
order - and a fortiori this is verified when all anisotropics take 
the value D' < D. Therefore Br-only DTN realizes a random 
arrangement of local environments (with D anisotropy or D' 
anisotropy) which are both gapless and ordered in the thermo- 
dynamic limit. 
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FIG. 3: Finite-size scaling of the correlation length at zero magnetic 
field and at T = 68.8 mK in the Br-only DTN model with J'^ = 
2.35 Jc and D' = 2Jc everywhere. 



The gapless nature of the spectrum of Br-only DTN in zero 
magnetic field is inherited by Br-doped DTN, which remains 
gapless down to zero magnetic field regardless of the doping 
concentration. This is due to the fact that Br-rich regions, lo- 
cally mimicking Br-only DTN, can have arbitrarily large sizes 
and therefore arbitrarily small local gaps. The Br-rich regions, 
although exhibiting local correlations similar to Br-only DTN, 
are not correlated with each other, so that the system in zero 
field does not possess long-range order, resulting in a Mott 
glass state. 



Local gap model 

The local gap model (LGM) describes a statistical ensemble 
of Br-rich clusters, that are gapless in the thermodynamic limit 
even in zero field. Exact diagonalization shows that small 
one-dimensional Br-rich clusters, modeled by the anisotropic 
Heisenberg model 

i i 

(with D' I J'c — 0.86) have a spectrum with an 7715- = ground 
state, and first excited states with ms = ±1 at an energy dif- 
ference A AT (where is the size of the cluster), vanishing 
in the thermodynamic limit. Therefore the application of a 
magnetic field to the system leads to a splitting of the first two 
excited states, leaving the ground state unaffected. The energy 
of the two excited states become then 

A^^ = (A^ T h) 0{An T h) (6) 

where h = gjiBH/ Jc. Here A^r is expressed in units of Jc- 
The 6>— function term is introduced ad hoc so as to make the 
cluster gapless when the local gap is closed by the field, and 
it prevents a gap from re-opening when h > A^. In fact a 
re-opening of the gap is not physically motivated, given that 
an increase of the field will make the size of the locally mag- 
netized region increase, thereby increasing the local density 
of states at low energy and not the opposite. For a cluster of 
size A^, we take the following scaling ansatz for the gap in a 
finite applied field 



which is typical of Heisenberg antiferromagnets undergoing 
spontaneous symmetry breaking [5]. Nonetheless, for irregu- 
larly shaped clusters this is obviously a strongly simplifying 
assumption, given that the above Ansatz only depends on the 
size of the cluster, and not on its geometry. 

At low temperatures, t = UbT/ Jc ^ ^(a^) with (A^) !^ 
1 + 1/xo, where xq — \ log(2x)| we can restrict the spectrum 
to the three low-energy states depicted above. 

Therefore the thermodynamics of the system reduces to that 
of a three-level system, with the following formulas for the 
thermal energy AEn{T) = En{T) - En{0) and total mag- 
netization Mn{T) of a cluster of size A^: 

M^(T) = ^ [e-z^^^^' - e-''^" '1 (8) 
Zn L J 

AEj,{T) = J-[AWe-'5<'+Aj,-'e-'^^^"'] (9) 
Zj, = l + e-^<+'+e-''<"' (10) 

Disorder averaging 

We now consider a generic observable On, that only de- 
pends on the size A^ of the cluster to which it is referred. The 
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average of the corresponding intensive quantity o over the en- 
semble of Br-rich clusters can be generally written as 



{o) 
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(11) 



N=l 



where ^at is the number of clusters of size N present in the 
system, can be expressed as the probability of having a 
given cluster of size at a given location in space and with 
a given orientation and geometry, times the number of differ- 
ent locations, orientations and geometries that the cluster can 
take. Defining a cluster of size N and delimited by M sites, 
we can then write 



9n = ^^^QnmPnm] 



(12) 



M 



here 



Pnm ~ exp(-xoA^) exp(-xiM) 



with xi = \ log(l — 2x)|; qnm represents the number of pos- 
sible geometrical realizations, rotations and translations of a 
cluster of size N and delimited by M sites. The translations 
are easily taken into account, by writing g^u = oinmL^- 
Therefore, we can generally write 



where 



gN = An exp{-xoN)L^ 



An = ^aArMexp(-xiM) 

M 



(14) 



(15) 



(we have reabsorbed the normalization of Pnm into the coef- 
ficient aNAi)- It is clear that the exact knowledge of An rep- 
resents a challenging problem. Therefore we make the strong 
simplifying assumption that An = A independent of N, and 
we use A as a fitting parameter. (Though we assume A is inde- 
pendent of A^, it still depends on the doping concentration x.) 
Our local-gap model depends therefore on two phenomeno- 
logical parameters only: c and A. The above assumptions 
lead to the expression 



(o);^A^O;vexp(-xoiV) (16) 

N=l 

which can be cast in an integral approximation 

(o) ^Af% 0{y) eM-^olv) (17) 

7o y 

justified by the fact that at low temperatures and low fields 
only large sizes contribute significantly to the thermodynam- 
ics of the system. The above integral can be easily computed 
numerically. 



Energy and specific heat in zero field 

We focus here on the thermodynamic properties in zero 
field, a limit in which the ensemble averages of the local gap 
model can be calculated almost exactly. The thermal energy 
(in units of Jc) of a cluster of size N in zero field reads 

^EN{t) = {c/N)/{l + exp(/3c/Ar)/2) (18) 

and the corresponding average thermal energy per spin reads 

dy e-^o/^ 



Ae{t) ^ 2cA [ 
Jo 



y 2 + 6^^2^ 



(19) 



A stationary phase approximation provides a very accurate ex- 
pression for the above integral in the following form 

2 \ 1/4 



Ae(t) ^2cA{^] t^^^ e^p(-2Jcxo/t) . (20) 

\xoc; 



(^^^ Consequently the specific heat reads 



C{T)/kB = 

A / 2 \ 1/4 

T (^) (l + 4V^)exp(-2V^) (21) 

giving Eq. (2) in the main text in the limit t -> 0. 



Comparison with QMC results 

Here we compare energy and magnetization obtained from 
the LGM model to the QMC results. The comparison is made 
between the thermal energy AE{T), and the thermal magne- 
tization Am(T) = m(T) — m(0). This is more appropriate 
than comparing the simple energy and magnetization. Indeed 
the LGM cannot predict the value of the ground state energy, 
and it only assumes the form a scaling Ansatz on the local 
gap. Moreover, the T = magnetization is not well captured 
by QMC at low fields, because it comes from extremely rare 
clusters whose size might be comparable with that of the sim- 
ulation box, and whose probability to appear, Pn, is too small 
to be correctly sampled by the disorder averaging performed 
in the QMC simulation. In fact, at T = the LGM model 
predicts the following expression for the magnetization: 



m{h]T = 0) 



2x0 



exp(— cxo/Zi). 



(22) 



This is the typical magnetization curve of the Bose glass [6, 7], 
which is exponentially suppressed at low fields. The essential 
singularity exhibited by the magnetization at h — implies 
that the Mott glass phase, realized in zero field, has a van- 
ishing linear susceptibility/compressibility, as well as vanish- 
ing susceptibilities to all orders of non-linearity, despite being 
gapless. 
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FIG. 4: Energy from LGM and QMC at various fields. 



We compare the QMC data with the predictions of the LGM 
by first fitting the parameters c and A to the QMC data for the 
energy at /i = 0. For the QMC parameters used in the simu- 
lation {x = 7.5%, Jc = 2.2 K, J ah = 0.18 K = 2.35 Jc, 
D — 8.9 K, D' — D /2), we find the best agreement between 
the LGM and the QMC data for c = 3.75 and A = 4.55. The 
/i = data for the energy are very well described by the LGM 
prediction. The values for A and c deviate from the ones that 
best fit the experimental results (c = 3.02 and A = 2.69) due 
to the strong sensitivity of the low-field thermodynamics to 
the actual value of the doping and of the Hamiltonian param- 
eters. This sensitivity comes from the fact that the dominant 
thermodynamic response of the system in low fields is given 
by the Br-rich regions, as captured by the local-gap model, 
and that the probability distribution of the size of such regions 
is exponentially dependent on the doping value. 

Remarkably, without any other adjustable parameter, we 
can extend the comparison between QMC and LGM to finite 
fields, finding an extremely good agreement up to fields h 1 
(corresponding to 1.5 T). The only significant deviations 
appear for h > hd, due to the appearance of the finite-T 
transition, which the LGM cannot account for. Moreover, the 
magnetization is also very well reproduced up to /i !^ 0.6, 
again without any further adjustable parameter. We attribute 
the systematic deviation between the low-temperature QMC 
data and the LGM data to the fact that the low values of 
the magnetization (< 10~^) are poorly sampled in the QMC 
simulations, due to the rarity of the regions with which they 
are associated, and to the finite size of the simulation boxes, 
which imposes an artificial upper cutoff on the size of the rare 
regions. 



The fact that the LGM model captures the thermodynamics 
of the Hamiltonian model for Br-doped DTN so well makes it 
a very flexible tool for the comparison with the experiments. 
The microscopic details of the Hamiltonian are absorbed in 




FIG. 5: Magnetization from QMC and LGM. 

the two fitting parameters c and A, while the doping and ap- 
plied field are external parameters that can be varied continu- 
ously. 



ESTIMATE OF THE CRITICAL FIELDS AND c 
EXPONENTS 
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FIG. 6: Windowing estimate of the exponent and of the critical 
fields from experiments (AC susceptibility measurements) and QMC 
simulations. 



The power-law scaling of the critical temperature as a func- 
tion of the critical field, - | - i/^ | ^ and Te - \H - 
Hc2\^, holds only in the vicinity of the critical fields Hd 
and Hc2, and it can be considered as the asymptotic limit of 
the field-dependence of the critical temperatures in the low- 
temperature regime. Its validity is evident in the log-log plots 
in Fig. (4) of the main text. To extract Hd, Hc2 and in the 
asymptotic limit one can adopt a windowing technique [8, 9], 



6 



which consists of approaching the low-temperature regime by 
fitting the experimental or numerical data over a temperature 
window T G [^.w] with a progressively smaller width. In 
practice we fit the temperature dependence of the critical field 
for the finite-temperature BEC, using the three parameter form 
H{T) = Hci(2) + ArV^i(2) in the vicinity of the critical 
fields, for temperatures in the [0, w] window. The results for 
both the exponents and critical fields are shown in Fig. 6. We 
see that both experiments and theory point towards a well- 
defined asymptotic regime with convergence of the critical 
fields and of the exponents towards a well-defined asymp- 
totic value, reported in the main text. The value and error bar 
of the fitting parameters are chosen as those associated with 
the w value at which the onset of convergence is observed. 
This allows us to exploit the largest number of data reproduc- 
ing the supposed asymptotic limit. 
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